/************************************************************************
int network_inferred (int N, int L, int pathway[L][N], int mode, int mode1, int mode2)

Input

N		number of nodes
L		number of steps
pathway		process
mode		0 : Process with steady state ; 1 : Cycle Process ; 2 : Triggled Cycle Process
mode1		0 : Budding Yeast Model ; 1 : Fission Yeast Model (with self-activation)
mode2		0 : Majority Rule; 1 : Restriction I ; 2 : Restriction II ; 3 : Restriciton I and II

Return Value	1 : Feasible

Note		Threshold Rule
***********************************************************************/

int Rij, Rii, Gij, Gii, theta;

int calc_design2 (int N, int L, int mode, int mode1, int mode2, int state[L][2 * N + 1], int i, int solution[N], int j)
{
	int k, l, NS, sum, flag;

	if (j == N)
	{
		for (k = 0; k < L; k++)
		{
			sum = 0;

			if (mode1 == 2 && k == 0)
				flag = 0;
			else if (mode2 == 0)
				flag = 0;
			else if (mode2 == 1)
			{
				if (state[k][N + i] == state[k][i])
					flag = 0;
				else
					flag = 1;
			}
			else if (mode2 == 2)
			{
				if (state[k][N + i] != state[k][2 * N])
					flag = 0;
				else
					flag = 1;
			}
			else if (mode2 == 3)
			{
				if (state[k][N + i] == state[k][i] && state[k][N + i] != state[k][2 * N])
					flag = 0;
				else
					flag = 1;
			}

			if (flag == 0)
			{
				for (l = 0; l < N; l++)
				{
					if (i == l)
						sum += solution[l];
					else
						sum += state[k][l] * solution[l];
				}
			}
			else
			{
				for (l = 0; l < N; l++)
				{
					if (solution[l] == 0)
						continue;
					else if (i == l)
						sum += solution[l];
					else if ((solution[l] > 0 && state[k][N + i] == 1) || (solution[l] < 0 && state[k][N + i] == 0))
					{
						sum += state[k][l] * state[k][N + l] * solution[l];
					}
					else if (state[k][l] + state[k][N + l] != 0)
					{
						sum += solution[l];
					}
				}
			}

			if ((sum > theta && state[k][N + i] == 1) || (sum < theta && state[k][N + i] == 0) || (sum == theta && state[k][N + i] == state[k][i])) 
				continue;
			else
				return 0;
		}

		return 1;
	}

	NS = 0;
	if (i != j)
	{
		solution[j] = Rij;
		NS += calc_design2 (N, L, mode, mode1, mode2, state, i, solution, j + 1);
		if (NS > 0)
			return 1;
		solution[j] = 0;
		NS += calc_design2 (N, L, mode, mode1, mode2, state, i, solution, j + 1);
		if (NS > 0)
			return 1;
		solution[j] = Gij;
		NS += calc_design2 (N, L, mode, mode1, mode2, state, i, solution, j + 1);
		if (NS > 0)
			return 1;
	}
	else
	{
		solution[j] = Rii;
		NS += calc_design2 (N, L, mode, mode1, mode2, state, i, solution, j + 1);
		if (NS > 0)
			return 1;
		solution[j] = 0;
		NS += calc_design2 (N, L, mode, mode1, mode2, state, i, solution, j + 1);
		if (NS > 0)
			return 1;
		if (mode1 == 1)
		{
			solution[j] = Gii;
			NS += calc_design2 (N, L, mode, mode1, mode2, state, i, solution, j + 1);
			if (NS > 0)
				return 1;
		}
	}

	return NS;
}

int calc_design (int N, int L, int mode, int mode1, int mode2, int pathway[L][N], int i)
{
	int j, k, state[L - 1][2 * N + 1], solution[N];

	for (j = 0; j < L - 1; j++)
	{
		for (k = 0; k < N; k++)
		{
			state[j][k] = pathway[j][k];
			state[j][N + k] = pathway[j + 1][k];
		}

		if (j != L - 2)
			state[j][2 * N] = pathway[j + 2][i];
		else
		{
			if (mode == 1)
				state[j][2 * N] = pathway[1][i];
			else
				state[j][2 * N] = pathway[j + 1][i];
		}
	}

	if (mode == 2)
	{
		if (i == 0)
			return 1;
		else
		{
			if (pathway[0][i] == pathway[1][i])
				solution[0] = 0;
			else if (pathway[1][i] == 1)
				solution[0] = Gij;
			else if (pathway[1][i] == 0)
				solution[0] = Rij;

			return calc_design2 (N, L - 1, mode, mode1, mode2, state, i, solution, 1);
		}
	}
	else
	{
		return calc_design2 (N, L - 1, mode, mode1, mode2, state, i, solution, 0);
	}
}

int network_inferred (int N, int L, int pathway[L][N], int mode, int mode1, int mode2)
{
	int i, NS, parameter;

	parameter = 3;
	if (parameter == 1) /* strong inhibition rule */
	{
		Rij = -2 * N;
		Rii = -1;
		Gij = 2;
		Gii = 1;
		theta = 0;
	}
	else if (parameter == 2) /* majority rule */
	{
		Rij = -2;
		Rii = -1;
		Gij = 2;
		Gii = 1;
		theta = 0;
	}
	else if (parameter == 3)
	{
		Rij = -2;
		Rii = -2;
		Gij = 2;
		Gii = 2;
		theta = 0;
	}
	else if (parameter == 4)
	{
		Rij = -2;
		Rii = -2;
		Gij = 2;
		Gii = 2;
		theta = 1;
	}
	else if (parameter == 5)
	{
		Rij = -2;
		Rii = -2;
		Gij = 2;
		Gii = 2;
		theta = -1;
	}

	for (i = 0; i < N; i++)
	{
		NS = calc_design (N, L, mode, mode1, mode2, pathway, i);
		if (NS == 0)
			break;
	}

	if (NS == 0)
		return 0;
	else
		return 1;
}
